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ABSTRACT 

We include an energy term based on Dark Matter (DM) self-annihilation during the 
cooling and subsequent collapse of the metal-free gas, in halos hosting the formation 
of the first stars in the Universe. We have found that the feedback induced on the 
chemistry of the cloud does modify the properties of the gas throughout the collapse. 
However, the modifications are not dramatic, and the typical Jeans mass within the 
halo is conserved throughout the collapse, for all the DM parameters we have con- 
sidered. This result implies that the presence of Dark Matter annihilations does not 
substantially modify the Initial Mass Function of the First Stars, with respect to the 
standard case in which such additional energy term is not taken into account. We 
have also found that when the rate of energy produced by the DM annihilations and 
absorbed by the gas equals the chemical cooling (at densities yet far from the actual 
formation of a proto-stellar core) the structure does not halt its collapse, although 
that proceeds more slowly by a factor smaller than few per cent of the total collapse 
time. 
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1 INTRODUCTION 

In the currently favoured ACDM cosmological model, the 
bulk of the matter component is believed to be made of 
(so far) electromagnetically undetected particles, commonly 
dubbed Dark Matter (DM). Although the evidence for the 
existence of DM is compelling on different scales, yet its 
nature is unknown, and many particle models beyond the 
standard one have been proposed in the literature as DM 
candidates. We address the reader to a recent review of ob- 
servational evidence and particle candidates for DM (e.g. 
Bertone, Hooper & Silk 2005), and will concentrate in this 
paper on a particular class of candidates, i.e. Weakly In- 
teracting Massive Particles (WIMPs). Many WIMP DM 
models are stable (under the conservation of the suitable 
symmetry, for each model) and hence do never decay into 
standard model particles. However, in many of these very 
same models the WIMPs are Majorana particles, thus car- 
rying the remarkable property of being self-annihilating; the 
value of the self-annihilation cross section, arising naturally 
in many WIMP models, reproduces the dark matter relic 



abundance required by the ACDM cosmology, if the mass 
scale of WIMPs is within the GeV/TeV scale and they are 
to be thermally produced in the early Universe. We adopt 
this as a benchmark scenario for our paper, and will often 
refer to it as a "Vanilla WIMP" . 

The actual DM distribution in the local Universe is such 
that even in the densest regions (e.g. galactic nuclei and 
black hole surroundings) from which the annihilation sig- 
nal could be in principle detected, the energy released by 
WIMP DM annihilations (hereafter, DMAs) is only a neg- 
ligible fraction of the one associated with standard gas pro- 
cesses. This implies that, locally, DM affects the host system 
almost uniquely through its gravitational effects, perhaps 
with the only possible exception of peculiar locations, such 
as the central parsec of the Milky Way (Fairbairn, Scott 
& Edsjo 2008, Scott, Edsjo & Fairbairn 2009, Casanellas 
& Lopes 2009). The effects of annihilating (or decaying, a 
scenario we do not consider here) DM upon the evolution 
of the intergalactic medium (IGM) at high redshift have 
been thoroughly studied (e.g. Chen & Kamionkowski 2004; 
Padmanabhan & Finkbeiner 2005; Mapelli & Ferrara 2005; 
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Mapelli, Ferrara & Pierpaoli 2006; Furlanetto, Oh & Pier- 
paoli 2006; Zhang et al. 2006; Ripamonti, Mapelli & Ferrara 
2007a; Valdes et al. 2007; Shchekinov & Vasiliev 2007), and 
are now believed to be small, except perhaps in the case of 
an extremely high clumping factor (Myers & Nusser 2008; 
Chuzhoy 2008; Natarajan & Schwarz 2008, Lattanzi & Silk 
2009), if one takes into consideration a standard, Vanilla 
WIMP scenario. 

The effects of DMAs upon primordial star-formation 
might be more significant. As the IGM could be heated by 
the energy deposition from DMAs, its temperature might in 
principle exceed the virial temperature of the smallest ha- 
los with the result of quenching gas accretion onto them. 
This effect has been shown to be unimportant by Ripa- 
monti, Mapelli & Ferrara (2007b, hereafter RMF07). How- 
ever, DMAs are expected to become more important as the 
collapse proceeds to protostellar scales (Ascasibar 2007). 
Spolyar, Freese & Gondolo (2008; hereafter SFG08) found 
that during the protostellar collapse of the first (Pop III) 
stars, the energy released by DMAs and absorbed by the 
gas could compensate (or even overcome) the radiative cool- 
ing of the gas. The increasing importance of such process 
arises from the combined enhancement during the collapse 
of DM density (due to gravitational dragging) and gas op- 
tical depth, implying a higher annihilation luminosity and 
absorption by the gas. The final phases of the collapse, after 
the formation of a hydrostatic core for gas central densities 
n c = pc/m v > 10 18 cm~ 3 (where p c is the central bary- 
onic density, and m p is the proton mass), have been inves- 
tigated by Iocco et al. (2008, hereafter 108), Freese et al. 
(2008b, 2009) and Spolyar et al. (2009). Initially, the DM 
pile-up is purely driven by gravitational interactions, but as 
the protostar approaches the Zero Age Main Sequence, DM 
accretion becomes dominated by the capture of WIMPs lo- 
cated in the star host halo after they scatter stellar baryons. 
As a consequence of the peculiar formation process of Pop 
III, following the smooth collapse of the gas cloud at the 
very center of the DM halos hosting them, Iocco (2008) and 
Freese, Spolyar & Aguirre (2008) suggested that DM capture 
is relevant for primordial stars; however, it can be safely ne- 
glected once local star formation is concerned, as the latter 
takes place anywhere in galactic discs, and it doesn't follow 
from a single, centered gas collapse episode. Further studies 
(108, Yoon, Iocco & Akiyama 2008; Taoso et al. 2008) have 
concluded that WIMP DM capture's most remarkable effect 
is the possible increase of the stellar lifetime. 

Quite surprisingly, the early phases of the collapse have 
received so far less attention with respect to the more ad- 
vanced ones, i.e. after hydrostatic core formation. For exam- 
ple, it is still unclear if the energy injection following annihi- 
lations results in a nett heating or cooling of the gas. In fact, 
high energy photons and electrons heat the gas through ion- 
izations; however, this heat input could be overwhelmed by 
the increased production of cooling species (as for example 
molecular hydrogen) stimulated by the larger abundances 
of free electrons, thus resulting in a net gas cooling. This, 
among others, is one of the aspects of the collapse of first 
stars in presence of WIMP annihilation that we would like 
to address here. We plan to do so by a set of sophisticated 
numerical simulations including all the relevant chemical re- 
actions and cooling processes. A first attempt to model the 
effects of DMA energy input was presented in Ripamonti et 



al. (2009); this paper represents a substantial extension and 
improvement of that study. 

Throughout the paper we assume the following set of 
cosmological parameters: Qa = 0.76, f2 m = 0.24, fib = 
0.042, fi DM = fim - fi b = fiwiMP = 0.198, and h = 0.73. 



2 METHOD AND CODE 

We base our investigation on a 1-D spherically symmetric 
code described by Ripamonti et al. (2002; hereafter R02). 
The original code, which includes the treatment of gravita- 
tion, hydrodynamics, and especially the chemistry and cool- 
ing of primordial gas, was originally conceived for the study 
of the last phases of the collapse of a primordial protostar 
(see also Omukai & Nishi 1998); later, it was extended in 
several ways (see Ripamonti 2007, Ripamonti et al. 2007a, 
and RMF07). 

Our simulations are based on those described in 
RMF07; here we list their most important properties, es- 
pecially when they differ from RMF07: 

• A single typical halo with mass 10 6 Mq virializing at 
z = 20 (virial radius i? v ir — 5 x 10 20 cm) is considered; the 
baryon fraction inside such halo is assumed to be equal to 
the cosmological value (fl^/QuM — 0.175); 

• The simulations are started at z = 1000 and involve a 
comoving volume 1000 times larger than that of simulated 
halo; initial baryonic density and temperature are constant 
in all the simulation shells, and equal to the cosmological 
values; 

• Before virialization, the gravitational effects of DM are 
treated as in the NFW case of RMF07: a predetermined (but 
time-dependent until virialization) DM potential is added to 
the gas self gravity. Such potential mimics the evolution of 
a halo in the top-hat approximation: as the DM potential 
becomes steeper, the (initially uniform) gas falls towards the 
centre of the halo, similarly to what is predicted by theory 
and consistent with the results of simulations (see e.g. Fig. 2 
of Abel et al. 2002); 

• After virialization, the evolution of the previously de- 
scribed artificial DM potential is stopped: its state at viri- 
alization is set to a NFW profile (Navarro, Frenk & White 
1996) with c = 10 and ifeoo = ^?vi r ; this is reasonably close 
to the results shown in Fig. 2 of Abel et al. (2002). The 
evolution of the DM profile is also followed in response to 
the baryon contraction (see later) in order to compute the 
DMA rat(Q 

• RMF07 investigated whether stellar formation might 
occur in a halo, whereas here we investigate how it starts; for 
this reason, simulations are stopped only when their compu- 
tational costs become very large (usually at number densities 
n c « 10 14 cm- 3 ); 

• The DM energy input is computed only after halo viri- 
alization, and only in regions with high baryon density 
[p > 4 x 10~ 22 gcm~ 3 , i.e. n = p/m p S 250cm~ 3 ), rather 
than at all times and everywhere: this is because RMF07 

1 We do not account for the gravitational effects of the adia- 
batically contracted DM profile, because adiabatic contraction is 
effective only when the baryonic potential largely dominates over 
the DM. 



already showed that before virialization and at low baryon 
densities the effects of DMAs are small; 

• For the purpose of evaluating the DMA rate (see be- 
low), the DM density Pdm(t~) is evaluated by assuming the 
conservation of the so called "adiabatic" invariant (see Blu- 
menthal et al. 1986). We implement the algorithm described 
by Gnedin et al. (2004), following the details in 108, and 
using the NFW and gas profiles described above as initial 
condition^; 

• The specific luminosity due to DMAs is 

Idm = <? (av) /)dm/(™dm), (1) 

where (av) is the thermally-averaged annihilation rate, and 
771dm is the WIMP mass; in the following we adopt (av) = 
3 x 10~ 26 cm 3 s _1 , whereas we consider 777dm as a free pa- 
rameter in the range 1.78 x 10" 24 g < m D M < 1.78 x 10~ 21 g 
(i.e. lGeV < TTiDMC 2 < ITeV); 

• The energy e that each baryon actually absorbs from 
DM annihilations (per unit time) is calculated through a 
detailed radiative transfer calculation, formally identical to 
the one performed for gray continuum radiation. Such calcu- 
lation is based on a constant gas opacity k — 0.01 cm 2 g" 1 
(roughly similar to the values used by SFG08). Moreover, 
since it is believed that the energy from annihilations splits 
roughly equally into electrons, photons, and neutrinos, we 
assume that only 2/3 of Ium (i.e. the fraction not going into 
neutrinos) can be absorbed; 

• Similar to RMF07, e can go into ionization, heating and 
excitation of atoms and molecules; we employ the results of 
Valdes & Ferrara (2008) (see also Shull & Van Steenberg 
1985, and Furlanetto & Johnson Stoever 2010) to estimate 
how to split the energy input into these three. Also note 
that the "ionization" component is split into ionization of 
H, D, He, He + , and dissociation of H2 , HD and Hj. In our 
"standard" treatment each species receives a fraction of the 
ionization energy which is proportional to its total baryonic 
content (in number, see RMF07 for details); 

• R02 switched to equilibrium chemistry (e.g. allowing 
the use of Saha equations instead of the detailed balance 
ones) for shells with number densities n > 10 13 ' 5 cm" 3 . Here 
we drop this simplification since (i) DMAs effects change 
the chemical evolution of the gas, and usually delay the ap- 
proach to equilibriurrjfl, and (ii) we never venture to densities 
n > 10 15 cm" 3 . 

Given the standing ignorance on the precise detail of 
feedback effects on the ionization and dissociation of atoms 
and molecules (especially H2), in addition to the standard, 
fiducial set of runs we performed runs with either enhanced 
or reduced feedback, in order to bracket the possible impact 
of such process. In the same way, since the opacity k we 



2 It is to be noticed that the DM profile in (gravitationally) 
baryon dominated regions is eventually dictated by the amount 
of gas accumulated after the contraction, see e.g. Fig 1 in SFG08. 
Our conclusions, especially regarding the final phases of the col- 
lapse, are hence almost insensitive to the initial DM profile. 

3 Even in the few cases where it is possible to switch to the 
equilibrium chemistry (e.g. the "control" run where we do not 
consider DMAs) we prefer to keep integrating the non-equilibrium 
equations in order to get results which are completely consistent 
with those from the other runs. 
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Table 1. Properties of the runs 



Name 


m D MC 2 [GeV] 


H 2 fbk 


Notes 


mi nnn 


1000 


std 


ivll Illlildl 


M100 


1 nn 


std 


Fiducial 


M10 


10 


std 


O LI VJ Hid Jvl 111 til 


Ml 




std 


1VI axi m al 




1000 


off 




N100 


i nn 


off 




N10 


10 


off 




Nl 


1 


off 




E1000 


1000 


enh. 




E100 


100 


enh. 




E10 


10 


enh. 




El 


1 


enh. 




L100 


100 


std. 


k = 0.001 cm 2 g" 1 


L10 


10 


std. 


k = 0.001 cm 2 g" 1 


LI 


1 


std. 


ft = 0.001 cm 2 g" 1 


H1000 


1000 


std. 


ft = 0.1 cm 2 g — 1 


H100 


100 


std. 


k = 0.1 cm 2 g — 1 


H10 


10 


std. 


ft = 0.1 cm 2 g~ 1 


NODM 






Control 



The first letter of the name of a run indicates the set to which it 
belongs; "M" refers to the main set, "N" to the set without DMA 
feedback upon H2 formation, "E" to the set with enhanced DMA 
feedback upon H2 formation, "L" to the set with low opacity, and 
"H" to the set with high opacity. The gray opacity is set to k=0.01 
cm 2 g — 1 for all of the runs presented in this Table, except for runs 
in the "L" and "H" sets. The NODM run assumes no energy input 
from DMAs. 

employ represents only a very rough estimate, we performed 
runs with either higher or lower values of k. 



3 RESULTS 

We test the effects of DMAs varying different sets of pa- 
rameters: (i) the normalization of the DMA rate, which is 
regulated by the ratio (av) /moM', (ii) the feedback strength 
on chemistry; and (iii) the "gray" gas opacity k. We antic- 
ipate that the strength of feedback on chemistry has little 
impact on the overall results, and that the effects of a vari- 
ation in k are somewhat similar to a variation of the same 
amount in (av) /moM- We will discuss the dependence on 
these two latter parameters in Section f3. 31 

Here we start by introducing in more details the physics 
of first star formation in presence of DMAs for our fidu- 
cial models ("M" labeled runs, see Table [T] for details). It is 
worth anticipating that the effects of DM become more rel- 
evant (and they become efficient earlier during the collapse) 
for higher self-annihilation cross sections/lower DM masses 
(see eq. [TJ. In what follows, we always adopt (av)—3x 10 -26 
cm 3 /s, and vary the particle mass tubm- Given the degen- 
eracy in the DMA energy injection term, the results can be 
interpreted at fixed mass and correspondingly varying the 
self-annihilation rate. 

Figure [1] shows the evolution of the temperature T c of 
the innermost shell as a function of the density of the same 
central shell of the simulated objects (n c ), for the five most 
representative models, M100, M10, Ml, M1000, and NODM 
(see Tabled. 
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Figure 1. Evolution of the temperature T c of the central shell 
of our simulation (whose mass is ~ 2 X 10 — 4 Mq) as a function 
of its baryon number density n c , starting slightly before virial- 
ization. The thick solid line shows the evolution in the control 
run with no DM energy input (NODM). The thin solid (blue) 
corresponds to the fiducial model (run M100); the dashed (red) 
to the maximal model (Ml); the dot-dashed (green) to the sub- 
maximal model (M10); and the thin dotted (magenta) to the min- 
imal model (M1000). 

As expected, DMA effects, which can be quantified by 
the deviations from the thick solid line, are most promi- 
nent in runs with low mass particles itiomc 2 = 1 GeV (Ml 
run, corresponding to the maximal DMA energy injection 
rate), and become very limited in the case of DM parti- 
cles with high mass m-DMC 2 = 1 TeV (M1000 run, minimal 
DMA energy injection rate); we will comment extensively 
on such dependence later in this Section. In the follow- 
ing we will refer to the "fiducial" run (M100; where the 
choice of parameters is quite standard), to the "minimal" 
run (M1000; where the high value of thdm reduces the en- 
ergy injection from DMAs, providing a check on the level 
where its effects become negligible), to the "maximal" run 
(Ml; where the parameters are chosen in a way which maxi- 
mizes the energy injection from DMAs and their effects upon 
the gas evolution), and to the "sub- maximal" run, M10. The 
maximal model, Ml, is likely to represent a sort of upper 
limit on the effects of DMAs on the formation of primor- 
dial stars. In fact, the DM parameters corresponding to this 
run ((av}=3x 10 -26 cm 3 s _1 , ?tidmc 2 = 1 GeV) are currently 
severely disfavored by multimessenger constraints on DM, 
in particular by (astrophysical) parameter-independent ones 
(e.g. GalH et al. 2009; Cirelli, Iocco and Panci 2009). 

3.1 The indirect feedback phase 

Figure [2] compares the total energy input from DMAs into 
the gas (which needs to be splitted into ionizations, heating 
and excitations) to the cooling rate due to H2 . This com- 
parison largely oversimplifies the thermal balance (see next 
subsection); however, it makes clear that, even in the maxi- 
mal case, the heating from DMAs remains a relatively small 
fraction of the H2 cooling until n c becomes > 10 9 cm -3 . At 
earlier stages (i.e. lower densities) the effects of DMAs are 
more subtle than what could be simply inferred by the direct 
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Figure 2. Evolution of the H2 cooling rate and of the rate of 
DMA energy input in the central shell of our simulations, as a 
function of its number density. In each panel the H2 cooling rate 
(thin solid line) and the DMA energy input (dashed line) for 
one of our runs (from top to bottom: fiducial, sub-maximal, and 
maximal runs) are compared to the H2 cooling rate in the control 
NODM run (thick solid line). We note that the actual heating 
rate from DMAs is about 1/3 of what is shown here, since the 
total DMA energy input also includes energy actually going into 
excitations and ionizations (see text for details). 

heating, which would not be able to significantly change the 
temperature of the gas. In fact, it is quite remarkable that 
in the 10 3 cm -3 < n c < 10 11 cm -3 range, the additional en- 
ergy injection from DMAs results in a decrease of the central 
temperature (see also Fig. [1]). 

Such counter-intuitive temperature decrease is due to 
indirect effects of DMAs, and in particular to their feedback 
on the chemistry; this remarkable property clearly emerges 
upon inspection of the curves in Fig. [2] showing that H2 
cooling (thin solid line) is stronger when DMAs are included 
(at least for n c < 10 10 cm" 3 ). 

In fact, as it can be inferred from Fig. [3] DMA ioniza- 
tion effects keep the gas more ionized than in the NODM 
case (i.e. e~ abundance 10 -6 -10 -4 , rather than <C 10 -6 ). 
In turn, this relatively high ionization level favors the for- 
mation of H2 through the reaction chain 

H + e~ -)■ H~+ 7 (2) 
H-+H Ha + e", (3) 

which is the main formation route for H2 in moderately 
dense, dust-free gas, and where the electrons act as cata- 
lysts (see e.g. Galli & Palla 1998; Glover & Jappsen 2007; 
Stasielak, Biermann & Kusenko 2007). 
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Figure 3. Evolution of the H2 and e~ fractions in the central 
shell of our simulation as a function of its baryon number density. 
In each panel the H2 fraction (thin solid line) and e~ fraction 
(thin dashed line) from one of our runs (from top to bottom: 
fiducial, sub-maximal, maximal runs) are compared to the same 
quantities from the control NODM run (thick solid line for H2 
and thick dashed line for e~ ). 



As Ripamonti et al. (2009) remarked, the increase in the 
H2 fraction can amount to 2-3 orders of magnitude, which 
apparently might imply same order-of-magnitude differences 
in the cooling rate, and large differences also in the gas tem- 
perature. However, while the increase in the cooling rate is 
quite large (though smaller than what is naively expected 
from the difference in H2 abundance), the corresponding re- 
duction of the gas temperature is much smaller, of the order 
of 30 per cent (see FigfTJ). This is due to the very strong tem- 
perature dependence of H2 cooling, especially below ~ 500 K 
(the temperature corresponding to the transition from the 
fundamental state to the lowest rotationally excited level): 
in fact, in this regime a moderate reduction in temperature 
results in a much larger decrease in the H2 cooling rate, 
which in turn slows down the temperature decrease. 

As a result, since the Jeans mass sensitivity to the tem- 
perature is quite mild, 

/ T \ 3/2 / n \-V2 
M j( n,T).5OM (_) (_) , (4) 

DMAs might reduce the Jeans mass scale but only by a 
factor < 2 (see Fig. Q}. 




4 6 8 10 12 14 

Log n c [crrT 3 ] 

Figure 4. Evolution of the Jeans mass calculated using the 
density and temperature in the central shell of our simulation, 
as a function of its number density. Lines refer to the minimal 
(M1000; dotted), fiducial (M100; thin solid), sub-maximal (M10; 
dot-dashed), maximal (Ml; dashed) models, and to the control 
run (NODM; thick solid). The top panel shows the absolute value 
of Mj, whereas the bottom panel compares the various models 
by showing the ratio of Afj in a given model to Mj ^jodm i n the 
NODM run. 



Table 2. Density and temperatures where the total energy input 
from DMAs overcomes H2 cooling. 



m DM c 2 


n crit [cm 


" 3 1 


T 

cr 


t[K] 


1 GeV( a > 


1.5 x 10 9 (~ 


l 9)(i>) 


780 (~ 


>800)<» 


10 GeV< a ) 


4.0 x 10 10 (~ 


l ii)W 


1070 r> 


1000) W 


100 GeV<» 


1.3 x 10 12 (~ 


io 13 )( fc ) 


1580 (- 


T300)( 6 ) 


1000 GeVM 


8.0 x 10 12 (~ 


1 W)W 


1740 (~ 


1400) W 



Notes: (a) The quantities in this table are almost independent 
of the assumptions about the strength of H2 feedback, so we do 
not distinguish among runs with the same moM- (b) Numbers in 
brackets indicate the results which can be estimated from Fig. 2 
of SFG08, intersecting the Yoshida et al. (2006) curve with the 
curves for various DM masses; for the SFG08 1 GeV mass case 
we adopt the 100% H2 line. 

3.2 The direct feedback phase 

Figure [2] shows that when the density increases above n c ~ 
10 12 cm -3 (~ 10 9 cm -3 ) for the fiducial (maximal) run, the 
energy input from DMAs finally overcomes the H2 cooling 
(see also Table [2] for details and a comparison with previous 
literature). This marks the beginning of what we call the 
direct feedback phase of the collapse, since the DMAs di- 
rect effects (especially the heating) finally start to dominate 
both H2 cooling and the more subtle DMAs indirect effects 
discussed above. We will refer to this condition (equality 
of DMA heating and H2 cooling terms) as to the "critical 
point" . 

The study of this phase is particularly interesting be- 
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cause when the DMA heating starts to compensate the ra- 
diative cooling, the protostar becomes unable to shed away 
its gravitational energy; it had been previously proposed (see 
SFG08) that this could induce the stop of the collapse, and 
the formation of a new type of stable astrophysical object 
powered by DMAs. 

The evolution of all the quantities beyond the critical 
point, shown in Figures [1] and [2] imply the existence of a 
dynamical evolution of the system also after the onset of the 
direct feedback phase -and hence no stalling of the object- 
even if there are important differences among the various 
runs. 

In the fiducial run, the evolution of the central shell 
proceeds relatively undisturbed. For n c < 10 14 gcm -3 the 
central temperature remains below the one of the control 
run, even if the gap tends to close (see Fig. [TJ . Conversely, 
in the maximal run there is a substantial steepening of the 
evolution in the n c — T c plane, and the central temperature 
overcomes the value of the control run at n c ~ 10 11 cm -3 . 
The temperature increase becomes particularly dramatic, 
with a sudden rise by a factor ~ 2 when the central density 
reaches n c ~ 10 12 cm -3 . 

In order to understand this phase, we need to describe 
the heating and cooling processes in detail. Fig. [5] shows 
the contributions from the main heating and cooling mech- 
anisms during the contraction of the central shell, for our 
reference runs. 



3.2.1 Details of chemical cooling 

Figure [5] shows the large importance of chemical heat- 
ing/cooling, i.e. of the energy which is released/absorbed 
as a consequence of chemical reactions, and in particular 
(in this regime, at least) of the formation/dissociation of 
H2 molecule^]: in every run there is a phase where H2 for- 
mation is the main heating mechanism. In the control run, 
chemical heating dominates over the adiabatic heating in 
the 10 9 ' 5 cm -3 H n c < 10 12 ' 5 cm -3 density range, which co- 
incides with the regime where 3-body reactions (see Palla, 
Salpeter & Stahler 1983) turn most of the Hydrogen into 
molecular form. In the runs with DMAs such a phase is 
even more extended, as the DMA feedback on the chem- 
istry anticipates the epoch of rapid H2 formation: in the 
maximal run, chemical heating becomes dominant already 
at n c > 10 3,5 cm -3 . 

It is remarkable that the chemical term can contribute 
both to heating and to cooling. It was already known (see 
e.g. Omukai & Nishi 1998; R02) that because of this prop- 
erty, it tends to act as a "thermostat" which stabilizes the 
evolution of the collapse, even if it cannot always prevent 
sudden transitions. 

The evolution of the maximal and sub-maximal models 
are particularly significant in this respect. In these models 
the H2 fraction peaks (at a level 0.8; see Fig. |3} when 

4 The formation of a H2 molecule releases its binding energy 
(4.48 eV); at the densities we are considering here, most of this 
energy is eventually converted into thermal energy of the gas (see 
e.g. Hollcnbach &; McKee 1979); on the other hand, the dissocia- 
tion of H2 through collisional reactions absorbs the same amount 
of thermal energy from the gas. 
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Figure 5. Evolution of the main heating and cooling mechanisms 
in the central shell of our simulations, as a function of its number 
density n c . The panels refer to four different runs: from top to 
bottom, we show results for the control (without DMAs), fidu- 
cial, sub-maximal, maximal runs. Each panel shows H2 line cool- 
ing (thin solid line), continuum cooling (thick solid line), DMA 
heating (dashed line), adiabatic heating (dotted line), and chem- 
ical heating or cooling (dot-dashed line). We note that i) the 
continuum cooling line can be plotted only when T c > 1000 K 
(corresponding to n c 10 9 — 10 1() cm - 3 ), however it is negligi- 
ble at lower temperatures; ii) the DMA heating is different from 
the DMA energy input shown in Fig. [2] since it does not include 
the ionization and excitation components; iii) the chemical term 
does not have a constant sign: shaded regions indicate the density 
ranges where it corresponds to a cooling term, whereas the ranges 
where it corresponds to a heating term are not shaded. 
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n c ~ 10 11 -10 12 cm' 3 : in such conditions 3-body H2 forma- 
tion is slowed down by the relative scarcity of atomic H left, 
and must compete with the increasingly large number of 
dissociations due to the ionization component of the DMA 
energy inpu10. 

When the molecular fraction starts to go down, the cool- 
ing rate must follow, because it is dominated by H2 and the 
number of molecules is being reduced; but so does the heat- 
ing, because its dominant component (chemical heating) is 
reduced by the energy lost in molecular dissociations. The 
fact that both the heating and the cooling mechanisms fol- 
low the same trend prevents a significant steepening of the 
temperature evolution and the consequent stop in the pro- 
tostellar collapse. This is because the "excess" energy which 
is not radiated away is mainly used to dissociate H2 . 

However, this cannot last indefinitely, as the reservoir 
of H2 is finite, and the DMA energy injection increases with 
increasing density. In the case of the maximal run, where 
the DMA energy injection is quite high, H2 is rapidly ex- 
hausted. At n c ~ 10 12 cm -3 there is essentially no H2 left, 
which also implies almost no radiative cooling. Then, the 
temperature goes up very fast (see Fig. [T]), developing a 
pressure gradient which is able to momentarily stop the col- 
lapse of the innermost shell (actually, there is even a very 
brief re-expansion episode- more details are discussed in Sec- 
tion [3]4]); however, the increase in temperature also results 
in a very strong increase in the continuum radiation which 
is produced by collision-induced emission, H~, and atomic 
H: continuum radiation increases by some 6 orders of mag- 
nitude (cfr. Lenzuni, Chernoff & Salpeter 1991; Marigo & 
Aringer 2009) and replaces H2 line radiation as the main 
cooling mechanism for the gas. The collapse is then able to 
resume, even if further investigations are needed in order to 
understand its details beyond n c » 10 13 cm -3 . At present 
in fact, computational costs become almost unaffordable be- 
yond n c ~ 10 13 cm~ 3 because of the strong sensitivity of 
cooling rate on temperature. In spite of these limitations 
we consider the conclusions based on these runs very robust 
as we have pushed the evolution to densities exceeding 10 4 
times the "critical" density at which the H2 cooling becomes 
comparable to the DMA energy transferred to baryons. 

Things are quite different in the case of the sub-maximal 
run: here a moderate H2 dissociation rate can compensate 
both the heating and the ionization parts of the (relatively 
small) DMA energy input for a significant span. In fact, 
at n c = 10 14 cm -3 (i.e. at a density ~ 100 times larger 
than the density at the onset of H2 dissociation, and ~ 1000 
times larger than the density at the "critical point" ) the H2 
fraction is still 0.1: it is still possible that H2 exhaustion 
(which we expect to happen at n c ~ 10 15 cm -3 ) is accompa- 
nied by a transition similar to the one we discussed for the 
maximal case; however, such transition is likely to involve 

5 It is important to note that the DMA energy injection which 
goes into ionizations actually contributes to the chemical heating. 
For example, if the number density of H2 molecules remains con- 
stant because the rate of formation (through 3-body reactions) is 
compensated by the rate of dissociations induced by the DMAs, 
the chemical heating term will be positive. In fact, the energy 
released by H2 formation ends up in the form of thermal energy, 
whereas the energy for the dissociations is not taken away from 
the gas thermal energy, as it is provided by DMAs. 
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Figure 6. Effects of H2 feedback schemes (upper panel) and opac- 
ity (lower panel): in the upper panel the solid (black) represents 
the maximal run Ml, the dashed (red) line is the El run (case 
with enhanced H2 feedback) and dotted (blue) line is the Nl run 
(no H2 feedback). In the lower panel thick (red) lines represent 
the Ml (solid) and Ll(dotted) runs, the intermediate (green) lines 
represent the the M10 (solid) the L10 (dotted) and H10 (dashed), 
finally, the thin (blue) line is our fiducial M100 run. See text for 
details and Table [T] for run coding. 

a smaller adjustment in temperature (say, an increase by a 
factor ~ 1.5 rather than ~ 2) because the required increase 
in the continuum cooling is much smaller (somewhat less 
than 3 orders of magnitude rather than 6, as can be inferred 
from Fig. [5)l than what was required in the maximal case. 

The fiducial and minimal runs exhibit trends which are 
somewhat similar to the sub-maximal run: the peak in the 
H2 abundance at n c ~ 10 12 cm -3 is followed by an extremely 
slow decline, during which the adiabatic heating is often as 
important as the DMA heating: then, it is quite likely that 
the transition from H2 line cooling to continuum cooling is as 
smooth as in the scenario without DMAs. In fact, it should 
be noted that at n c = 10 14 cm -3 the continuum cooling 
amounts already to about 10 per cent of the total heating, 
and it is increasing faster than any heating mechanism (see 
Fig. El). 

3.3 Sensitivity to feedback and opacity 

In order to check the effect of our assumptions about H2 
feedback (in particular, destruction of H2 molecules by the 
energy injected by DMAs into the gas) we have run two extra 
sets of simulations: the first (label: E) in which H2 feedback 
is enhanced, and the second (label: N) with DM-induced H2 



s 



dissociations switched off. In the E-models, we have assumed 
that all DMA energy input which goes into excitations is in 
the form of photons in the Lyman- Werner band (11.2 eV < 
hv < 13.6 eV), and that each of these photons dissociates 
one H2 molecul^fl This clearly represents an upper limit on 
the feedback upon H2 , as in reality about half of the photons 
going into dissociations has energy below the lower limit of 
the Lyman- Werner band; furthermore, in low density regions 
a significant fraction of Lyman- Werner photons can escape. 

As it can easily be noticed from the upper panel of Fig- 
ure [6l (showing the evolution of models Ml, El and LI in 
the n c — T c plane) switching off the H2 feedback results in 
a sensible change of the temperature only in the final stages 
of collapse (because the absence of DMA-induced H2 disso- 
ciations strongly reduces the heating); on the other hand, 
strengthening the feedback results in appreciable modifica- 
tion of the temperature of the cloud in earlier stages (be- 
cause it reduces the H2 abundance), but the difference re- 
duces drastically in later phases (when the direct effects of 
DMAs dominate the ones on chemistry). We also note that 
in models with ttidmc 2 > 10 GeV the differences induced 
by the treatment of feedback are much smaller than in the 
models shown in Fig. [6] (with jtidmc 2 = lGeV). 

We have also studied the effects of varying the gray 
opacity k, by running sets of simulations in which we have in- 
creased/decreased the fiducial value k=0.01 cm 2 g _1 by one 
order of magnitude. At low and intermediate densities the ef- 
fects on the evolution of the cloud are very similar to those of 
a decrease/increase of the particle mass didm; whereas mod- 
els with the same ihdm tend to converge at high densities. 
In fact, an opacity variation affects sensibly the properties 
of the gas only when it is of "optically thin" to the energy 
produced by DMAs: in later stages of the collapse, (i.e. in 
regions of the halo where the gas is "optically thick" to the 
DMA energy), the energy absorption does not depend on 
k. This can be appreciated from the lower panel of Figure 
[6] where it is clear that the modified opacity models (H, L 
labeled) overplot to the corresponding standard (M labeled) 
models at higher central densities. 

We can definitely infer that "astrophysical" (and "nu- 
merical") parameter dependence is mostly degenerate with 
DM parameters, and in any case not drastically affecting the 
qualitative picture we have drawn so far. 

3.4 The evolution of spatial profiles 

3.4-1 The fiducial run 

Figure [7] shows the spatial profile of gas properties at four 
different stages (corresponding to central shell gas densities 
n c = 10 5 , 10 8 , 10 11 , and 10 14 cm -3 ) intheMlOO fiducial case 
and, for comparison, in the NODM control run. The top-left 
panel shows also the DM density in the case of M100 run 

6 In practice, we assume a Lyman- Werner H2 dissociation rate 
(per unit volume) equal to ne/ exc / -Elw 1 where n is the baryon 
number density, e is the energy input per baryon due to DMAs, 
/exc is the fraction of this energy which goes into excitations (i.e. 
into photons with 10.2 eV < hu < 13.6 eV), and i?LW — 12 eV is 
the average energy of a Lyman- Werner photon. Such rate is added 
to the one produced by the fraction the DMA energy input which 
goes into ionizations. 
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Figure 7. Evolution of DM and gas profile (top left), infall 
velocity (top right), temperature (bottom left), and H2 fraction 
(bottom right) during the collapse, comparing runs M100 (fidu- 
cial) and NODM (control). Solid lines refer to gas quantities in 
the M100 case, dashed lines to gas quantities in the NODM con- 
trol run. Different solid lines refer to central shell gas densities 
n c =10 5 (black), 10 s (blue), 10 11 (green), 10 14 cm" 3 (red); the 
different stages can be easily distinguished through the position 
of their innermost point, as the lines referring to high-density 
stages extend to lower radii. In the upper left panel, the dotted 
lines refer to DM densities in the M100 case (note that the DM 
density profile for the n c = 10 5 cm -3 stage is very similar to the 
NFW profile which is assumed to form immediately after virial- 
ization) . The geometrical markers indicate the position enclosing 
baryonic masses of 10 4 Mq (hexagons), 1O 2 M0 (pentagons), IMq 
(squares), 1O _1 M (triangles), 10~ 2 M Q (circles), and 10~ 3 M© 
(diamonds); see text for comments. 



at the same stages (DM density stays constant and equal to 
the initial profile in the NODM control run). 

The first important remark to file is that the gas density 
profiles are virtually unchanged by the presence of DMAs: 
such profiles are indistinguishable until central gas densities 
of n c ~10 8 cm -3 are reached, and even then the discrep- 
ancies are of minor entity (e.g., the slope outside the core 
becomes ~ 2.15 rather than ~ 2.20). We can conclude that 
DMAs do not alter the self-similar phase of the collapse 
(Larson 1969, Penston 1969), at least for run M100 and for 
n c < 10 14 cm~ 3 . The evolution of the DM density profile 
is in good agreement with previous studies, see for instance 
Fig. 1 in SFG08, and Fig. 1 in 108. 

It is also relevant to notice that (as could be expected) 
the differences in infall velocity, temperature profile and 
H2 fraction affect only the central ~ 10 4 Mq of gas (i.e. 
the gas within the hexagonal marks). Actually, relevant (al- 
though yet not dramatic) differences are limited to the inner 
~ 10 2 Mq (region within the pentagonal marks): the overall 
properties of the 10 6 Mq halo are conserved in presence of 
DMAs. 

The lower panels visualize the changes in temperature 
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Table 3. Collapse times (yr) of the innermost gas shell 



Log (n c [cm 3 ]) 


NODM 


M100 


Ml 


2-3 


2.117X10 6 


2.116x10 s 


2.107x10 s 


3-4 


8.382x10 s 


8.326x10 s 


8.077x10 s 


4-5 


3.204x10 s 


3.200 x10 s 


3.401x10 s 


5-6 


1.158x10 s 


1.177x10 s 


1.190x10 s 


6-7 


3.930 xlO 4 


4.026X10 4 


3.895X10 4 


7-8 


1.261 xlO 4 


1.250X10 4 


1.221X10 4 


8-9 


3.927X10 3 


3.757X10 3 


3.544 xlO 3 


9-10 


1.210x10 s 


1.118X10 3 


1.072x10 s 


10-11 


366.4 


333.3 


333.7 


11-12 


103.2 


102.0 


111.3 


12-13 


29.42 


31.93 


31.05 


13-14 


9.210 


9.480 





The time intervals needed for the innermost shell to increment 
its density of one order of magnitude, in the NODM, fiducial 
and maximal run. The highest density is absent in the Ml case 
because the run has been stopped at smaller density, see text. 

(lower left) following the modification of H2 fraction (lower 
right) which we have extensively commented in the previous 
sections. Here we only note that the molecular core is sig- 
nificantly more extended than in the control run: the region 
where H2 is the dominant chemical species encloses ~ 5 Mq , 
rather than ~ 1 Mq ; furthermore, the fall in H2 abundance 
is less steep than in the case without DMAs. 

The upper right panel is somewhat intriguing, as it gives 
a glimpse of the modification of the dynamical properties 
of the collapse: while there are very little changes in the 
central region (say, within the square marks, corresponding 
to a 1 Mq enclosed mass) and in the outskirts of the halo 
(at enclosed masses jj 10 3 M©, i.e. outside a point roughly 
mid- way between the pentagonal and hexagonal marks), the 
infall velocity of the intermediate-mass gas shells is slowed 
down. The reduced temperature of the cloud (due to the 
more efficient H2 cooling which we discussed in the previous 
sections; see also the bottom-right panel of Fig. [7J) is likely 
responsible for this change. In fact, it is well known (see 
e.g. Stahler, Palla & Salpeter 1986, and references therein) 
that during the self-similar evolution, the mass infall rate 
pO")t>infali(r) in the region where the density falls as a power 
law with index ~ 2 should roughly scale with c 3 /G, where 
c a oc T 1//2 is the isothermal sound speed. Since p(r) is practi- 
cally the same in the M100 and the NODM run, the reduced 
temperature should result in a reduced infall velocity, similar 
to what we obtain. 

It is worth noting that on the contrary, the col- 
lapse of the innermost gas shell is accelerated, albeit very 
slightly, by the presence of DM, until central densities of 
10 12 cm~ 3 (/10 u cirr 3 ) for the M100(/M1) case, see Tabled] 
However, such infall time alteration are of very small entity, 
and in any case negligible with respect to the total collapse 
time - of the order of 3 x 10 s yr. 



3.4-2 The maximal run 

Figure [S] analogue to Fig. \7\ compares the control run to 
the maximal run, rather than to the fiducial run; in this case 
the highest central density we plot is n c = 10 13 cm -3 (rather 
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Figure 8. Same as in Fig. [7]but comparing runs Ml (maximal) 
and NODM (control). Solid lines refer to gas quantities in the Ml 
run, dashed lines to gas quantities in the NODM control run. All 
the symbols and lines are the same as in Fig. [7] apart from the 
fact that the red lines refer to a central gas density n c =10 13 cm -3 
(rather than 10 14 cm - 3 ); see text for comments. 



than 10 14 cm -3 ), because the maximal run was stopped be- 
fore reaching n c = 10 14 cm -3 . 

The top-left panel shows that even in this case DMAs 
do not alter the self-similar part of the collapse; however, it 
should be remarked that the highest density stage is starting 
to deviate from the expected self-similar profile: in fact, the 
density profile within the core is not flat, but there is a 
maximum at the edge of the core. 

This is more clear when we examine the other panels: 
the lines referring to the three stages with n c < 10 11 cm -3 
are qualitatively similar to their counterparts in Fig. [7] but 
the lines representing the highest density stage are clearly 
different, and we will focus on them. 

In the lower panels we can see a huge drop in the H2 
abundance near the centre of the protostar, and a corre- 
sponding discontinuity in the gas temperature; both could 
be expected by our discussion in the previous sections, even 
if here we get a further piece of information about the size 
of the region where the H2 was dissociated (R < 3 x 10 _s pc, 
enclosed mass ~3x 10 -2 Mq). Also the velocity panel is af- 
fected: this is the only case where the infall velocities of the 
run with DMAs significantly exceed those in the control run 
in an extended region; furthermore, we can see that there 
are secondary peaks in velocity at low radii (in particular 
one is coincident with the edge of the region where H2 was 
dissociated), which are likely to result from the propagation 
of the density wave which was generated by the re-expansion 
of the core which occurred when n c ~ 10 12 cm -3 . 

All these properties indicate that the core of this object 
is approaching hydrostatical equilibrium, even if it is not 
quite there. It is interesting to notice that the innermost 
regions, of enclosed mass between IMq and IOMq, are in- 
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falling with higher velocity than in the NODM case, whereas 
between lOM© and 1O 3 M0 the collapse is slowed, when cen- 
tral density is n c = 10 12 cm" 3 . Again, this is equivalent to 
what we have observed for the fiducial run, and shown in 
Table [3j and the feature only becomes more pronounced in 
the case of the maximal run. 



4 DISCUSSION 

We have studied the effects of self-annihilating Dark Matter 
on the collapse of the gas structures harboring the forma- 
tion of the first stars, known as Pop III. For the first time 
in the literature we follow self-consistently the evolution of 
the Dark Matter profile as a consequence of the gravita- 
tional drag of the collapsing gas and include the feedback 
of energy injection by DM annihilations on the chemistry 
of the gas, in the yet unexplored regime between the virial- 
ization of the halo and the formation of a hydrostatic core. 
We have explored the effects of DMAs by spanning a range 
of masses and annihilation cross sections around the values 
of the Vanilla WIMP scenario, namely those able to repro- 
duce the relic abundance of DM with a thermal decoupling, 
finding similar results but with variations in the details and 
onset times of the different phases that we are to summarize. 

In the following, quoted numbers refer to our fidu- 
cial case (m D MC 2 = 100 GeV, (av) = 3 x 10~ 28 cm 3 s _1 , 
k = 0.01 cm 2 g _1 ). As expected, heavier WIMPs (or smaller 
self-annihilation rates, or lower opacity since the energy in- 
jection term is proportional to (av) /dtdm and depends also 
on k) effects become relevant at later times during the col- 
lapse, the same DMA rate being achieved at higher DM 
densities. Lighter particles (or higher (av), or higher k) pro- 
duce the same effects at earlier times (and therefore smaller 
gas densities). This general scenario is robust with respect to 
parameter variations within physically acceptable ranges. A 
few key points of our study are worth some emphasis: (i) As 
known from existing literature (RMF07), before the virial- 
ization of the halo, at gas density n c < 10 3 cm -3 (and there- 
fore extremely small gas opacity) DMAs do not sensibly af- 
fect any gas process; (ii) between halo virialization and a gas 
density n c ~ 10 11 cm -3 DMAs contribute mainly through 
indirect feedback effects: the free electron floor created by 
the ionizations induced by DMAs catalyzes H2 formation. 
In turn, molecular H2 provides more cooling to the cloud 
than in the standard case (without DMAs) and the temper- 
ature of the cloud decreases as a consequence of DMAs; (Hi) 
finally, at n c > 10 11 cm -3 the DMA heating rate becomes 
equal to the gas cooling rate. To a first approximation this 
leads to a balance between losses and gain, which we dubbed 
as the "critical point". 

Our results generally agree with previous ones: as the 
semi-analytical estimates of SFG08 and the analytical study 
(based on simulation data) from Natarajan, Tan and 'O Shea 
(2009), we confirm the existence of a "critical point". We 
also find that the equality of H2 cooling vs DMA heating 
terms takes place (in the innermost shell) at central den- 
sities of approximately n c ~ 10 12 cm" 3 , details depending 
on parameters (DM mass, opacity coefficient etc. see previ- 
ous section, and following discussion), in agreement with the 
above mentioned studies. This is particularly relevant as our 
analysis is the first fully numerical, self-consistent study in- 



cluding a reasonably accurate treatment of radiative transfer 
allowing to reach the critical point and beyond it in at least 
three cases. 

It had been previously suggested (SFG08) that, after 
reaching of the "critical point" , the collapse would stop and 
the whole structure stall, thus generating a new type of celes- 
tial object. With our numerical simulations we have accessed 
this stage of the collapse, finding that the system does not 
halt its collapse in three cases out of four; in the only case 
where the collapse stops, this happens far after the critical 
point was reached, and its duration is very short (~ 3yr) 
after which the gas restarts contracting. Moreover, the DM 
parameters for which the astrophysical system finds impor- 
tant changes after the critical point are actually strongly 
disfavored by DM constraints based on local and primordial 
Universe observations. By changing the DM mass or self- 
annihilation rate, the scenario we have described does not 
change qualitatively, within the range of values studied and 
the physical regime we have accessed with our simulations. 

While stressing again that when our objects reach the 
critical point the structure does not halt the collapse, and 
instead it continues its evolution toward the formation of a 
protostellar core, it is also useful to comment on the alter- 
ation of the dynamics of the collapse induced by DMAs. 

The duration of the gas collapse from halo virialization 
down to densities of n c ~ 10 14 cm -3 , in presence of DMAs 
appears to be shorter (especially in the first phases) than 
in the standard case (by about 1 per cent) of the total col- 
lapse time. However we point out that the change is a small 
fraction of the total, and that the effects of DMAs on col- 
lapse time are very difficult to quantify. In fact, the DMAs 
effect can indeed accelerate the collapse in some shells, and 
decelerate it in others, thus making it difficult to obtain a 
homogeneous definition of time delay. For example, the re- 
duction of the infall velocity of shells enclosing ~ 10 2 Mg 
implies that the collapse time of these shells is longer (by 
10-20 per cent) than in the control (NODM) case. 

We feel confident in stressing, however, that the shells 
which are mostly affected by the infall time decrease (~20% 
of the infall velocity) are placed between 1O 2 M0 and 1O 3 M0, 
and would become part of a hydrostatic protostar only in 
the final stages of the collapse. It is not clear whether these 
shells will eventually end up in the hydrostatic core (Omukai 
& Palla 2003), however it is reasonable to expect that such 
modification will not alter the total time for the formation 
of the star by more than ~20 per cent. 

4.1 Hints about the hydrostatic core 

Even if our runs could not reach the regime when a hydro- 
static core forms (with the partial exception of the maximal 
run, where we probably got close), it is worth to examine the 
likely consequences of our results for the further evolution 
of the protostar in presence of DMAs. 

We start reminding that in the standard case (without 
DMAs) the simulations of R02 show that when the pro- 
tostar becomes opaque to continuum radiation (at n c ~ 
10 16 cm -3 ), the formation of the hydrostatic core is delayed 
by the thermostatic effects of H2 dissociations, so that the 
core actually forms when n c 10 20 cm -3 . 

However, at the end of our runs the central H2 abun- 
dance in the maximal and sub-maximal runs is < 0.1; even 
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in the fiducial run, the H2 abundance, while still high, is 
decreasing much earlier than in the standard case; this is 
very different from the results of R02, where /h2 ~ 1 up 
to n c > 10 16 cm -3 . Because of the lower amount of H2, it 
is reasonable to expect that in the cases with DMAs the 
delaying effects of H2 dissociation are absent, or smaller 
than in the standard case. Then, the hydrostatic core would 
form at lower densities (e.g., n c ~ 10 16 cm~ 3 , if the proto- 
star becomes optically thick to continuum radiation at the 
same density as in the R02 runs), and its initial mass would 
be larger (probably in the 0.01-0.1 Mq range, rather than 
~ 0.003 M s ). 

However, we point out that R02 found that in the case 
with no DMAs the mass of the hydrostatic core grows very 
fast, reaching 0.1 Mq in less than 3 years (see e.g. their Fig. 
6): then, the difference in the initial size is relatively unim- 
portant. It is probably more relevant to note that the lower 
temperatures and infall velocities of the layers outside the 
hydrostatic core imply that in the cases with DMAs the ac- 
cretion rate might be slower than what was found by R02. 



5 CONCLUSIONS 

We have studied the effects of WIMP Dark Matter Anni- 
hilations (DMAs) on the evolution of primordial gas clouds 
hosting the first stars. We have followed the collapse of gas 
and DM within a 1O 6 M0 halo virializing at redshift z = 20, 
from z = 1000 to slightly before the formation of a hy- 
drostatic core, properly including gas heating/cooling and 
chemistry processes induced by DMAs, and exploring the 
dependency of the results on different parameters (DM par- 
ticle mass, self-annihilation cross section, gas opacity, feed- 
back strength). Independently of such parameters, when the 
central baryon density, n c , is lower than the critical density, 
ro C rit ~ 10 9_13 cm -3 , corresponding to a model-dependent 
balance between DMA energy input and gas cooling rate, 
DMA ionizations catalyze an increase in the H2 abundance 
by a factor ~ 100. The increased cooling moderately reduces 
the temperature (by « 30 per cent) but does not signifi- 
cantly reduce the fragmentation mass scale. For n c > n cr it, 
the DMA energy injection exceeds the cooling, with the ex- 
cess heat mainly going into H2 dissociations. In the pres- 
ence of DMA the transition to the continuum dominated 
cooling regime occurs earlier and generally is not associated 
with abrupt temperature variations. In conclusion, no signif- 
icant differences are found with respect to the case without 
DMAs; in particular, and contrary to previous claims, the 
collapse does not stall and the cloud keeps contracting even 
when n c 3> Merit • Our simulations stop at central densities 
~ 10 14 cm -3 , and cannot follow the hydrostatic core for- 
mation, nor its accretion. At the final simulation stage, the 
lower temperature/infall velocity of the layers enclosing a 
mass of ~ 10 2 Mq suggest that DMAs might lead to slightly 
longer stellar formation timescales, with a possible « 20 per 
cent increase over models without DMAs. The latter find- 
ing strengthens our conclusions, although the final answer 
will come from numerical simulations (hopefully also three- 
dimensional) able to address this very same problem in the 
yet unexplored density regime 10 14 cm -3 < n c ~ 10 18 cm -3 . 
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APPENDIX A: APPLICABILITY OF THE 
ADIABATIC CONTRACTION FORMALISM 

One of the hypothesis underlying the AC formalism (as de- 
rived by Blumenthal et al. 1986) is that the orbital time of 
DM particles should be short when compared to the collapse 
time of the baryons. In the case we are studying, this is a 
problem: in fact, the baryonic collapse happens essentially 
on a free-fall time-scale, and the DM orbital times cannot 
be shorter than that. 

The case of a "fast collapse" (where the baryon collapse 
time is similar to the free-fall time) was studied by Steigman 
et al. (1978; see their sec. Illb), in the assumptions that i) 
the baryons stand still at t = (so that they will collapse 
in free-fall), ii) the DM contribution to the total mass M is 
small, and iii) DM particles have a dispersion velocity Av^, M . 
The case can be treated analytically, and Steigman et al. 
(1978) obtain that if the radius of the sphere enclosing all the 
baryons goes from Ri (at time t = 0) to R c (at time t = t c ), 
the ratio of the number JVdm, c of DM particles within R c 
(at t c ) to the number A^dm,i of DM particles within Ri (at 
t = 0) is roughly proportional to (R c /Ri) 3 ^ 2 (cfr. eqs. 31 
and 32 of Steigman et al. 19783), at least in the limit of 
Rc <C Ri- Then, the density pdm of DM particles should 
grow asymptotically as (R c / Ri)~ 3 ^ 2 , or as pj/ 2 . 

We can compare this behaviour to the dependence of 
the central DM density (pdm,c) upon the central baryonic 
density (pb,c) which can be obtained from the AC formalism: 
SPG08 find that pdm, c oc p^ 1 , and our calculations essen- 
tially confirm their result. Then, it is clear that the "fast 
collapse" slope is significantly natter than what can be in- 
ferred from the application ot the AC formalism. However, 
there are some extra facts that need to be kept into account. 

• For central baryonic number densities n c < 10 5 cm -3 
the DM density is larger than the baryonic one, and the ap- 
plication of the AC does not modify the DM density profile; 



7 Note that he power law index in eq. 32 of Steigman et al. (1978) 
should be —3/2 rather than —1/2, since this is the value that can 
be inferred from the more general eq. 31. 



• The "fast collapse" slope was obtained with the as- 
sumption of an uniform density profile of both the baryons 
and the DM, whereas the NFW profile of our initial 
DM profile is strongly peaked: then, the calculations of 
Nr>M,c/Nr>M,i by Steigman et al. (1978) likely underesti- 
mate the value appropriate in our case; as a consequence, 
it is likely that the slope of the pdm — Pb relation is higher 
than 0.5; 

• The DM profile seen in the simulations of Abel et al. 
2002 appears to be more consistent with the AC predic- 
tions than with those of the "fast collapse" model: in fact, 
the innermost DM density shown in their Fig. 2 is a fac- 
tor ~ 10 lower than the baryonic density n ~ 10 9 cm -3 : 
this is somewhat consistent with the difference by a factor 
[n/(10 5 cm" 3 )] 1-0 ' 81 ~ 6 predicted by AC, but quite smaller 
than the difference by a factor [n/(10 5 cm" 3 )] 1 " ' 81 ~ 100 
predicted using the "fast collapse" . 

Then, it appears that the AC formalism can be safely 
used at least for densities up to n c ~ 10 9 cm" 3 . It is quite 
possible (and even likely) that at higher density it will sig- 
nificantly overestimate the DM density. However, such over- 
estimate will be large (a maximum < 30 for n c = 10 14 cm -3 ) 
but still reasonable in the range of densities where we per- 
form our simulations. Furthermore, we remark that the AC 
approach allows a better comparison with previous studies, 
and, given our conclusions, is the conservative choice. 



